Segregation of morphogenetic regulatory function of Shox2 from its cell fate guardian role in sinoatrial node development

Shox2 plays a vital role in the morphogenesis and physiological function of the sinoatrial node (SAN), the primary cardiac pacemaker, manifested by the formation of a hypoplastic SAN and failed differentiation of pacemaker cells in Shox2 mutants. Shox2 and Nkx2-5 are co-expressed in the developing SAN and regulate the fate of the pacemaker cells through a Shox2-Nkx2-5 antagonistic mechanism. Here we show that simultaneous inactivation of Nkx2-5 in the SAN of Shox2 mutants (dKO) rescued the pacemaking cell fate but not the hypoplastic defects, indicating uncoupling of SAN cell fate determination and morphogenesis. Single-cell RNA-seq revealed that the presumptive SAN cells of Shox2-/- mutants failed to activate pacemaking program but remained in a progenitor state preceding working myocardium, while both wildtype and dKO SAN cells displayed normal pacemaking cell fate with similar cellular state. Shox2 thus acts as a safeguard but not a determinant to ensure the pacemaking cell fate through the Shox2-Nkx2-5 antagonistic mechanism, which is segregated from its morphogenetic regulatory function in SAN development.

In contrast, Tbx18 deficient mice form a hypoplastic yet functional SAN 12 .These threads of evidence implicate that the molecular mechanisms underlying SAN pacemaking cell fate determination and morphogenesis, a biological process that governs the development of tissue or organ shape by orchestrating the spatial arrangement of cells during embryonic development 21 , are not identical.
Previous studies showed that Shox2 is indispensable for the physiological function and morphogenesis of the SAN, as Shox2 null mice display hypoplasia of the venous pole structures and failed differentiation of SAN cells, embryonic bradycardia, and upregulation of chamberprogram genes in the venous pole structures [22][23][24] .Nkx2-5 was originally thought to be detrimental to SAN development, and Shox2 was shown to regulate SAN development by preventing ectopic Nkx2-5 activation in the SAN head domain 10,22,23,25,26 .Further studies showed that Shox2 counterbalances Nkx2-5 in the developing SAN junction domain, an interface between the SAN and atrium that is specialized by Nkx2-5 expression, to regulate the cell fate decisions between pacemaking cells and working myocardial cells 27 .We have recently found that Nkx2-5 inactivation in the SAN junction cells confers them with the property of the SAN head cells and impairs SAN function, but does not affect SAN morphogenesis 28 , further supporting the existence of morphogenic genetic programs in the SAN that are separated from pacemaking cell fate determination.However, the molecular mechanisms that distinctly control each process remain unknown.

Results and discussion
Maintenance of pacemaking cell fate in the developing SAN requires a minimal level of Shox2 expression To better understand how Shox2 controls SAN morphogenesis and differentiation, we investigated the physiological consequence of reduced Shox2 expression at different dosages during SAN development.We used a Shox2-Cre knock-in allele that express Cre recombinase instead of Shox2 as Shox2 null allele in our studies 29 .The Shox2-Cre heterozygous mice appeared indistinguishable from wildtype mice, never exhibited any obvious unusual behaviors, abnormal ECG recordings, and altered gene expression profiles.In the Shox2 HA-Neo/+ allele that was used as a hypomorphic allele in this study, the presence of the neomycin cassette interferes with Shox2 expression, resulting in a reduction of Shox2 expression to 65% 30 (Supplementary Fig. 1a; Supplementary Data 1).A combination of different Shox2 alleles enabled the manipulation of Shox2 expression at different dosages (Supplementary Table 1).As Hcn4 reflects the pacemaking capability of cells in the SAN, we examined Hcn4 expression in mice carrying different compound Shox2 alleles at various stages.We found that Shox2 HA-Neo/Cre mice, which had Shox2 production reduced to 15% and ectopic Nkx2-5 expression (Fig. 1b) compared to wildtype mice (Fig. 1a), displayed Hcn4 + domains in the Shox2 + lineage derived venous pole structures (Fig. 1b).However, such mice had a reduced SAN size with retained head and junction domains and shortened venous valves (Fig. 1b).These results indicate that a 15% level of Shox2 production compared to wildtype was sufficient for the determination and maintenance of pacemaking cell fate, but such residual dosage failed to support normal SAN morphogenesis.As compared to its requirement for normal SAN morphogenesis, the maintenance of pacemaking cell fate requires a lower dose of Shox2.The hypoplastic SAN morphology but with the retaining of the pacemaking cell fate in the Shox2 HA-Neo/Cre suggests the existence of distinct Shox2-mediated genetic modules in SAN development dictated by different Shox2 levels.
Deletion of Nkx2-5 in Shox2 mutant background rescues the pacemaking cell fate and ameliorates the impaired physiological function of the SAN We have previously reported that Shox2 maintains the pacemaker program in the Nkx2-5 + junction domain of the SAN by antagonizing Nkx2-5 transcriptional output 27 , we therefore asked whether rebalancing the Shox2-Nkx2-5 antagonism can reset the original pacemaking cell fate in the developing SAN.We took advantage of the compound Shox2 Cre/Cre ;Nkx2-5 F/F alleles to delete Nkx2-5 in the Shox2 null background 29,31 .In contrast to the embryonic lethality of Shox2 Cre/Cre mice at mid-gestation due to cardiac defects, Shox2 Cre/Cre ;Nkx2-5 F/F mice (dKO) survived to postnatal day 0 (P0), indicating amelioration of the impaired physiological function of the SAN.Surface electrocardiogram (ECG) measurements on P0 mice showed regular but slower cardiac rhythm and the normal P waves in the dKO mice as compared to controls (n = dKO: 7/7; WT: 16/16) (Fig. 2a).Additionally, whole-cell patch-clamp recording on cells derived from the Shox2 + population in the SAN of E13.5 dKO mice showed typical AP configurations comparable to controls (Fig. 2d; Supplementary Fig. 2; Supplementary Data 1).Immunohistochemistry and in situ hybridization assays showed that, unlike the lack of Tbx3 and Hcn4 expression in Shox2 null mice 22 , dKO mice retained the expression of SAN markers, including Tbx3, Hcn4 (Fig. 2e).This observation strongly implies the preservation of the pacemaking cell fate in dKO mice.Morphometric analyses on the SAN of dKO and Shox2 null mice (Supplementary Fig. 1b-d) revealed that the SAN morphology in dKO mice appeared comparable to that in Shox2 null mice 22,23 .However, upon closer inspection, we discovered that the Hcn4 expression in the SAN head domain (as opposed to the junction domain demarcated by Nkx2-5 expression) of dKO mice was notably restored (Fig. 2b; Supplementary Fig. 1b-d), following the deletion of Nkx2-5 in the Shox2-null background.This result suggests that the dKO mice retain the pacemaking cell fate by developing a hypoplastic SAN lacking the junction domain (Fig. 2c; Supplementary Data 1).Given previous reports indicating that the structural characteristics of the microenvironment surrounding pacemaking cells in the SAN can impact their function through local mechanics 32 , the sinus bradycardia ECG defects observed in dKO mice (Fig. 2a) may be attributed to the altered SAN morphology (hypoplasia).Above all, these observations demonstrate that rebalancing Shox2-Nkx2-5 antagonism by deletion of Nkx2-5 in Shox2-null background rescues the pacemaking cell fate and ameliorates the impaired SAN physiological function but has no impact on SAN morphogenesis observed in Shox2 -/- mice.The exclusive function of Shox2-Nkx2-5 antagonistic mechanism in the regulation of SAN cell fate, but not in SAN morphogenesis, suggests that there are alternative genetic programs utilized by Shox2 to control SAN morphogenesis.
and performed single-cell RNA-seq (scRNA-seq) on Shox2 + cells (labeled by GFP) from the SAN of E13.5 mice with these three different genotypes (control, Shox2-null, and dKO).The datasets obtained were subjected together to unsupervised clustering, resulting in 8 distinct clusters, including Fibroblast, endothelial cells, mesothelial cells (MCs), cardiomyocytes (CM), Macrophages, and epithelial cells (EPIs) 28,[33][34][35][36][37][38][39] (Fig. 3a).Based on the expression of Tnnt2 and GFP, we were able to identify a Shox2 + cardiomyocyte population (Fig. 3b, c).The Shox2 + cardiomyocytes were then further distinguished into three sub-populations (sP0, sP1, and sP2) (Fig. 3d; Supplementary Fig. 3a).While sP0 and sP2 showed the closest phylogenic distance and expressed pan-SAN marker Hcn4 40 (Fig. 3e, f), these two subpopulations were defined as pacemaker cells, the sP1 sub-population seemed likely the residual venous valve cells, as they have similar transcriptome with the venous valve we reported previously 28 .Further unsupervised sub-clustering reclassified the pacemaker cells into two groups, with the control and dKO pacemaker cells falling in group 0 (G0) and Shox2-null pacemaker cells in group 1 (G1) (Fig. 3g; Supplementary Fig. 3b, b').These results indicate that the control and dKO cells share a high transcriptomic similarity that deviates from that of Shox2-null cells.To infer the cellular states of the pacemaker cells, we used Monocle2 41 to pseudotemporally order the pacemaker cells along a developmental trajectory to the differentiated pacemaking cells.Pseudotime trajectory analysis showed that while the control and dKO samples exhibit the most differentiated cellular state, Shox2-null sample displays the most primitive and underdifferentiated cellular state (Fig. 3h, i; Supplementary Fig. 3c).This is evidenced by the increased expression of several widely recognized cardiac progenitor cell markers 42 , including Abcg2.Abcg2, which serves as the molecular marker for identifying side population progenitor cells in multiple adult tissues 43 , including the adult heart 44 , is also documented to be transiently expressed during the early stages of heart development 45 .Thus, the pacemaker cell fate-oriented differentiation process appears unaffected in the dKO cells but is blocked at a progenitor state in the Shox2-null cells, as compared to the control group.These results are in line with the rescued pacemaking cell fate in the dKO mice as shown above, providing a solid foundation for further analyses.
With the distinct physiological and morphological traits of the SAN in the three types of mice with distinct genetic background mentioned above, we reasoned that comparative analyses of the transcriptomic profiles between different sample would allow us to segregate the morphogenetic regulatory factors that are downstream of Shox2 from the pacemaking cell fate determination factors (Supplementary Table 2).To identify genes that are potentially involved in SAN morphogenesis, we compared the transcriptomic profiles between the dKO and control mice using Seurat 46,47 .We identified 1025 genes that are differentially expressed between the dKO and control mice (set 1).Gene ontology (GO) analysis on the set 1 differentially expressed genes (DEGs) revealed the enrichment of morphogenesisassociated GO terms, such as muscle tissue morphogenesis, muscle organ morphogenesis, and cardiac muscle tissue morphogenesis (Fig. 4a).In an attempt to identify genes that are involved in SAN cell fate determination, we then compared the transcriptomes between the dKO and Shox2 null samples and identified 1462 DEGs (set 2).Conversely, the identified set 2 DEGs included genes that are known to participate in pacemaking cell fate regulation, such as Hcn4, Smoc2, and Gata6 48 .These results suggest the feasibility of segregating the SAN morphogenetic regulatory factors from the pacemaking cell fate determination factors by comparative analyses on the transcriptomic profiles of the three samples.To segregate them, we then intersected the DEGs in set 1 and set 2 and found out that except 561 genes that were shared by these two sets, there are 464 genes that are exclusively presented in set 1, and 901 genes that were found exclusively in set 2 (Fig. 4c).These two gene sets are thus potentially involved in the segregation of SAN morphogenesis from cell fate determination.Heatmap analysis confirmed the differential expression of these genes (Fig. 4b).Furthermore, the differential expression of several selected genes including Hcn4, Smoc2, Bmp10, Gata6, Mef2a, and Tbx18 was further confirmed by pseudotime trajectory analysis and by in situ hybridization and immunostaining assays (Fig. 4d, e).Interestingly, although it was shown previously that Pitx2 and Tbx20 act upstream of Nkx2-5 22,30,49,50 , we observed ectopic expression of Pitx2 and Tbx20 in the scRNA-seq profile of the dKO SAN (Fig. 4b), suggesting that the Shox2-Nkx2-5 antagonistic mechanism has a negative feedback effect on the upstream regulator genes, which warrants future investigation.Above all, our unique genetically modified mouse lines allowed us to segregate genes involved in the regulation of morphogenesis from those involved in cell fate determination downstream from Shox2 and Nkx2-5 in SAN development.
ChIP-seq analysis reveals that the segregation of SAN morphogenesis from cell fate determination attributes to a functional interplay between Shox2 and Nkx2-5 Since the two gene sets are identified from Shox2 and Nkx2-5 mutant mice, to verify the identified genes and to further understand how Shox2 interacts with Nkx2-5 to regulate SAN morphogenesis and cell fate determination separately, we re-analyzed our published Shox2 and Nkx2-5 ChIP-seq datasets on the developing heart 27 .As outlined in our previous report 27 , Shox2 and Nkx2-5 exhibit genome-wide co-occupancy (Supplementary Fig. 4a).However, characterization of binding peak distribution patterns for Shox2 and Nkx2-5 unveils that Shox2 exhibits a higher affinity for binding to the promoter regions (82.65%) compared to Nkx2-5 (64.81%).Conversely, Nkx2-5 displays a stronger preference for binding to the distal cisregulatory elements, such as distal intergenic regions (Supplementary Fig. 4b).Intersection of Shox2 and Nkx2-5 binding peaks reveals that the cobinding regions of Shox2-Nkx2-5 are predominantly located within the promoter regions, similar to the distribution of Shox2 peaks (Supplementary Fig. 4c).Interestingly, when functionally categorizing the target genes associated with Shox2 binding, Shox2-Nkx2-5 co-binding, and Nkx2-5 binding peaks, it becomes evident that these genes play distinct roles in cardiac-related GO terms.Specifically, genes bound by Shox2, including those that are co-bound by Shox2-Nkx2-5, were observed to contribute to pacemaker cell differentiation and SAN morphogenesis (Fig. 5a).On the other hand, the genes exclusively bound by Nkx2-5 were linked to the cardiac cell fate commitment (Fig. 5a).This result aligns with Nkx2-5's binding preference and is consistent with our previous study, where we highlighted its role in defining a subpopulation of pacemaker cells 27 .More interestingly, only the genes with individual binding of Shox2 and Nkx2-5, rather than co-binding of Shox2-Nkx2-5, are classified within the context of cardiac pacemaker cell development and differentiation (Fig. 5a), implicating that the segregation of Shox2's morphogenetic regulatory function from its cell fate safeguarding role arises due to combinational interactions between Shox2 and Nkx2-5.These results not only reinforce our segregated gene outcomes, but also underscore the significance of delving deeper into the intricate mechanisms of Shox2-Nkx2-5 antagonism.Therefore, we conducted an intersection between the list of 3999 genes displaying exclusive Nkx2-5 binding (Fig. 5b) and the gene set 2 (Fig. 4c), which is potentially involved in cell fate determination.This analysis aimed to pinpoint genes associated with SAN cell fate commitment, yielding a list of 333 genes directly targeted by Nkx2-5.Notably, this list of genes encompasses those that have been previously implicated in cardiac cell fate determination, such as Smoc2 48 (Fig. 5c) and Taz 51 (Supplementary Data 2).To distinguish genes involved in safeguarding cell fate from those associated with morphogenesis, we performed further intersections on the following gene groups: firstly, the group of 1242 genes with both Shox2 and Nkx2-5 binding, but no Shox2-Nkx2-5 co-binding; and the group of 4691 genes with all three binding patterns of Shox2 and Nkx2-5 (Fig. 5b).These intersections were performed with gene set 1, yielding a list of 105 genes possessing both Shox2 and Nkx2-5 binding sites but lacking co-binding sites (Fig. 5d, Supplementary Data 3), and a separate list of 429 genes with all three binding patterns (Fig. 5e, Supplementary Data 4).While the 429 genes with all three binding patterns are potentially involved in SAN morphogenesis, we assumed that these 105 genes with both Shox2 and Nkx2-5 binding but no co-binding sites act as safeguarding factors.Interestingly, among these 105 genes, we found the DEGs that were shown in the regulation of both cell fate determination and morphogenesis through scRNA-seq analysis, including Tbx18, Pitx2 (Supplementary Fig. 4d, e), and Sox11(Fig.5d).This further substantiates the notion that the segregation during SAN development results from the functional interplay between Shox2 and Nkx2-5.Taken together, our integrative analyses of ChIP-seq and scRNA-seq datasets verified and identified genes that are distinctly involved in the regulation of pacemaking cell fate commitment, safeguarding and SAN morphogenesis, respectively.However, since the ChIP-seq assays were performed on the whole atrium of the heart that contains the SAN, atrial myocardial cells, and non-myocardial cells, while the scRNA-seq was done on the SAN cells, further studies are needed to verify our identified genes.
In this study, we have uncovered that the essential role of Shox2 in SAN morphogenesis is uncoupled from its function in promoting pacemaking cell fate.We have unraveled that the altered pacemaking cell fate in the SAN of Shox2 mutant mice can be reset by rebalancing Shox2-Nkx2-5 antagonism in Shox2 and Nkx2-5 dKO mice, despite the presence of a hypoplastic SAN.A morphologically integrated SAN is essential for mammals to retain the dominant pacemaker site in a proper position 52 .Theoretical studies and  empirical evidence point toward a morphogenetic program for the SAN that is evolved independently from the pacemaker program 12,53,54 .In this study, we showed that despite the presence of well-differentiated pacemaking cells in the SAN, the dKO mice nevertheless displayed bradycardia symptom, which was also seen in Shox2 null mutants, indicating the impact of the altered morphology on the physiological function of the SAN.Our scRNAseq on the SAN cells from those unique mouse lines allow us to establish a Shox2 and Nkx2-5 regulated gene expression profile exclusively implicated in the regulation of SAN morphogenesis along with a gene expression profile specifically involved in pacemaking cell determination and differentiation, which warrants future investigation and verification.These gene expression profiles would also provide candidate genes for the etiology of congenital and acquired abnormalities of the SAN and potential targets for therapeutic intervention of sinus node dysfunction.

Fig. 1 |
Fig. 1 | A minimal level of Shox2 expression maintains pacemaking cell fate in the developing SAN.Immunofluorescence shows the expression of Hcn4, Nkx2-5, and Shox2 in the SAN head, SAN junction, and venous valves in E14.5 wildtype control (a) and Shox2 HA-Neo/Cre mice (b) (n = 3 for each genotype).Note in the whole SAN (magenta dash-lines), SAN head (red arrowheads), SAN junction (green arrowheads), and venous valves (yellow arrows), compared to controls (a), Shox2 HA-Neo/Cre mice (b) exhibited residual Shox2 expression, but visible ectopic Nkx2-5 expression in the SAN head domain and reduced size of these structures.However, Hcn4 expression was retained.RA right atrium, RSVC right superior vena cava.Bar: 50 μm.

Fig. 3 |
Fig. 3 | ScRNA-seq identifies transcriptomic profile of dKO SAN cells comparable to controls.a UMAP visualization of the eight clusters from the SAN and its adjacent atrial cells of E13.5 Shox2 Cre/+ ;R26R mTmG (control), Shox2 Cre/Cre ;R26R mTmG (Shox2 null), and Shox2 Cre/Cre ;Nkx2-5 F/F ;R26R mTmG (dKO) mice.b, c VlnPlot shows the expression of Tnnt2 and GFP probability distributions across eight clusters.d-f Analysis of three sub-populations from Shox2 + (represented by GFP + ) cardiomyocytes defines sP0 and sP2 as the SAN cells.g Sankey Plot showing the SAN cells of control, Shox2 null, and dKO reclassified into two groups.h, i Pseudotemporal ordering of SAN cells reveals that Shox2 null cells exhibited the pseudo-start state, while the control and dKO cells displayed the pseudo-end state.Cells are colored by the genotype of samples (i).Cells with dark color represent the pseudo-start, and the bright color represent the pseudoend (h).CM cardiomyocyte, EC endothelial cells, MC mesothelial cells, EPI epithelial cells, C5, C7 cluster number 5, 7.

Fig. 4 |
Fig. 4 | Validation of selected transcription targets of Shox2 involved in SAN morphogenesis v.s.pacemaking cell fate determination.a Go analysis on the DEGs between control and dKO reveals GO terms associated with heart development.b Heatmap shows the genes involved in the determination of the pacemaking cell fate and SAN morphogenetic regulation across the three samples.c Venn diagram depicts the number of overlap genes between pacemaking cell fate determination and SAN morphogenetic regulation.d Pseudotime analysis unravels elevated expression of Gata6, Hcn4, Mef2a, and Smoc2 along pseudotime trajectory.e In situ hybridization and immunostaining verified the expression of Smoc2, Tbx18, Bmp10, and Hcn4 in the SAN of E13.5 control, Shox2 null, and dKO mice.RA right atrium, RSVC right superior vena cava.Bar: 50 μm.

Fig. 5 |
Fig. 5 | Integrative analyses identify direct targets of Shox2 and Nkx2-5 involved in the segregation of pacemaking cell fate commitment, safeguarding and SAN morphogenesis.a Functional categorization of Shox2 and Nkx2-5 binding genes.b Venn diagram depicts the number of overlapped genes between Shox2 binding, Nkx2-5 binding and Shox2-Nkx2-5 co-binding genes.c-e UCSC genome browser track displaying the binding patterns on the representative genes that are potentially involved in the segregation of pacemaking cell fate commitment, safeguarding, and SAN morphogenesis.